# Code for Table M1, Row 5: Permutation Test

rm(list = ls())

## ---------------------------------------
## Load Packages 
## ---------------------------------------
library('AER')
library('ivpack')
library('data.table')
library('doParallel')
library('parallel')

## ---------------------------------------
## Load Data and Functions
## ---------------------------------------
load("../data0812.RData")

directory <- "../functions/"
functions <- list.files(directory)  
loadfunctions <- sapply(functions, FUN = function(x)source(paste0(directory, x)))

## ---------------------------------------
# Construct Instrument
## ---------------------------------------
data0812 <- constructIV(data0812)

## ---------------------------------------
#  Select Last Case Before Election
## ---------------------------------------
data0812 <- lastCase(data0812)

## ---------------------------------------
## Two Stage Least Squares Specifications
## ---------------------------------------
time.controls <- "as.factor(court_time1) + as.factor(court_time2) + as.factor(court_dow) + as.factor(court_shift) + as.factor(totOGS2)"
case.controls <-   "as.factor(any_drug_2) +  as.factor(any_violent_2) + as.factor(fire_arms_2) +  as.factor(any_rob_2) + as.factor(any_dui_2) + as.factor(prior_offender_2)"
demo.controls <- "age_2012 + I(age_2012^2) + Female + as.factor(race) + vote2008 + as.factor(noteli08) + regis_before"

outc.1 <- "vote2012"
endo.1 <- "pti"
inst.1 <- "judgeiv"

form.1 <- formula(paste(outc.1, "~", endo.1, "+" , time.controls, "|", inst.1, "+", time.controls))
form.2 <- formula(paste(outc.1, "~", endo.1, "+" , time.controls, "+", demo.controls, "|", inst.1, "+", time.controls, "+", demo.controls))
form.3 <- formula(paste(outc.1, "~", endo.1, "+" , time.controls, "+", demo.controls, "+", case.controls, "|", inst.1, "+", time.controls, "+", demo.controls, "+", case.controls))

m1a1 <- ivreg(form.1, data = data0812)
m1a2 <- ivreg(form.2, data = data0812)
m1a3 <- ivreg(form.3, data = data0812)

cont.o <- data.table(round(rbind(c(robust.se(m1a1)[2, c(1:1)], robust.se(m1a2)[2, c(1:1)], robust.se(m1a3)[2, c(1:1)])
), 3))

cont.o

nc <- detectCores()
np <- 500
cl <- makeCluster(nc)
registerDoParallel(cl)

results <- foreach(i = 1:np) %dopar% {
  
  load("../data0812.RData")
  
  library('AER')
  library('ivpack')
  library('data.table')
  library('stringr')
  
  set.seed(i)
  assignments <- list()
  list.m <- 1:6
  t <- sample(list.m, 3) 
  
  for(m in 1:max(data0812$ashift, na.rm = T)) {
    assignments[[m]] <- t
    list.m <- list.m[-which(list.m %in% t)]
    if(length(list.m) == 0) {
      list.m <- 1:6
    }
    t <- sample(list.m, 3) 
  }
  
  data.assign <- data.table(do.call('rbind', assignments))
  data.assign$ashift <- 1:nrow(data.assign)
  
  data0812 <- merge(data0812, data.assign, by = c("ashift"))
  
  data0812$judge_cat <- NA
  data0812$judge_cat[data0812$court_shift == 1] <- data0812$V1[data0812$court_shift == 1]
  data0812$judge_cat[data0812$court_shift == 2] <- data0812$V2[data0812$court_shift == 2]
  data0812$judge_cat[data0812$court_shift == 3] <- data0812$V3[data0812$court_shift == 3]
  
  data0812 <- constructIV(data0812)
  
  data0812 <- lastCase(data0812)
  
  ## ---------------------------------------
  ## Two Stage Least Squares Specifications
  ## ---------------------------------------

  time.controls <- "as.factor(court_time1) + as.factor(court_time2) + as.factor(court_dow) + as.factor(court_shift) + as.factor(totOGS2)"
  case.controls <-   "as.factor(any_drug_2) +  as.factor(any_violent_2) + as.factor(fire_arms_2) +  as.factor(any_rob_2) + as.factor(any_dui_2) + as.factor(prior_offender_2)"
  demo.controls <- "age_2012 + I(age_2012^2) + Female + as.factor(race) + vote2008 + as.factor(noteli08) + regis_before"
  outc.1 <- "vote2012"
  endo.1 <- "pti"
  inst.1 <- "judgeiv"
  
  form.1 <- formula(paste(outc.1, "~", endo.1, "+" , time.controls, "|", inst.1, "+", time.controls))
  form.2 <- formula(paste(outc.1, "~", endo.1, "+" , time.controls, "+", demo.controls, "|", inst.1, "+", time.controls, "+", demo.controls))
  form.3 <- formula(paste(outc.1, "~", endo.1, "+" , time.controls, "+", demo.controls, "+", case.controls, "|", inst.1, "+", time.controls, "+", demo.controls, "+", case.controls))
  
  m1a1 <- ivreg(form.1, data = data0812)
  m1a2 <- ivreg(form.2, data = data0812)
  m1a3 <- ivreg(form.3, data = data0812)
  
  cont <- data.table(round(rbind(c(robust.se(m1a1)[2, c(1:1)], robust.se(m1a2)[2, c(1:1)], robust.se(m1a3)[2, c(1:1)])
  ), 3))
  cont
}

stopCluster(cl)

results.final <- do.call('rbind', results)

cat("\nPrinting Table M1, row 5: Permutation Test...\n")
out <- rbind(cont.o, cbind(mean(abs(results.final$V1) < abs(as.numeric(cont.o[1,1]))),
                    mean(abs(results.final$V2) < abs(as.numeric(cont.o[1,2]))),
                    mean(abs(results.final$V3) < abs(as.numeric(cont.o[1,3])))))

colnames(out) <- c("(1)", "(2)", "(3)")
out$reported_value <- c("Estimate", "Empirical p-value")
print(out[, c("reported_value", "(1)", "(2)", "(3)"), ])
